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ABSTRACT 



A mathematical model for the simulation of non-catalytic 
solid-gas reactions in a fluidized bed reactor is proposed. 
The performance of the model in predicting solid reactant 
conversions for an ore roasting process is investigated 
using available literature data. Model development required 
simplifying assumptions. The sensitivity of the model to 
certain of these assumptions is investigated. 

Comments on the adaptability of the model for use in 
the design and study of a fluidized bed shipboard waste dis- 
posal system are made. 
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I. INTRODUCTION 



In the late sixties, the Navy began a study of shipboard 
waste disposal systems which could be used to combat a grow- 
ing pollution problem. One proposal considered a process 
for the combustion of solid waste material in a fluidized 
bed reactor. This proposal was opposed by some naval 
designers on the grounds that sufficient knowledge of the 
operational characteristics of fluid bed reactors was not 
available, thus prohibiting accurate design and control of 
such a system. 

This lack of fundamental design knowledge is a problem 
which has plagued fluid bed technologists for some time. 
Although the fluid bed reactor finds broad usage in the 
chemical industry today, its complex nature makes accurate 
mathematical modeling difficult. Design and control of 
fluid bed systems in the past has, therefore, been based on 
the application of operational engineering experience. 

This work was undertaken with the goal of studying and 
developing a mathematical model of a fluidized bed reactor. 
The Bubble Assemblage Model proposed by Wen and Yoshida [1] 
provided the starting point for this investigation. 

If a reasonable model could be developed, it could be 
used as a design tool and in studies of the operational 
characteristics of the reactor system. 
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II. LITERATURE SEARCH 



Early tracer studies [2] on large scale fluid bed reac- 
tors revealed that tracer concentrations in the exit gas 
stream exceeded concentrations in samples taken directly 
from the bed. This observation led to the development of 
the two phase model of fluidization which predicts that the 
major portion of the gas passes through the bed in the form 
of bubbles. The. primary developers of this theory include 
Toomey and Johnstone [3], Shen and Johnstone [4], Pansing 
[5] , and Lewis et al [6] . 

The two phase model pictures the bed as consisting of 
two single phase reactors in parallel as shown in Figure 1. 
Phase 1, the lean or bubble phase, represents the primary 
means of gas throughput. It is generally treated as being 
in plug flow and containing no solids. Phase 2, the dense 
or emulsion phase, contains the solids surrounded by inter- 
stitial gas. Assumptions on the mixing patterns in the 
emulsion phase range depending on the investigator. Mathias 
and Watson [7], Massimillia and Johnstone [8], and Gomez- 
plata and Shuster [9] assumed plug flow behavior while Shen 
and Johnstone [4] and Lewis et al [6] have studied condi- 
tions of plug flow and complete mixing. The model is also 
characterized by a gas exchange coefficient which accounts 
for gas flow between the phases. 

May [10] extensively studied the classic two phase model 
and was the first to introduce diffusion type coefficients 
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Figure 1. Six Parameter Two Phase Model. 
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to characterize the vigorous mixing in the emulsion phase. 
This approach was investigated further by van Deemter [11] 
who applied the model in a study of gas mixing experiments 
in order to determine the diffusion and gas exchange 
coefficients . 

The mathematical formulation of the two phase model leads 
to two partial differential equations containing six para- 
meters and representing the system material balances. A 
general expression of the model equations can be given as 

p9CGB -nrvu ^ CGB 9CGB ^ mr-o rnT:'\ r> r i 

FDb ^^2 F^(CGB-CGE) = 0 (1) 



9CGE 



£De + £Ue + F (CGE-CGB) = 0 

9t2 9H o 



( 2 ) 



Their solution has been a center of controversy because of 
different approaches to the choice of boundary conditions. 
McCracken [12] obtained solutions for various mixing 
patterns and in a second work [13] reviews other numerical 
solutions of the model equations. 

The two phase model, although representing a useful 
engineering approach to the understanding of fluid bed 
reactors, suffers from serious disadvantages. The model 
parameters are, in a sense, adjustable and, therefore, while 
fitting experimental data well for particular cases, their 
generality for use in scale up and design are questionable. 
Furthermore, the description of a uniform bubble phase and 
constant gas exchange coefficient ignores effects of chang- 
ing bubble size on bed operation. 
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Pioneering work in development of a more general model 
was completed by Davidson and Harrison [14] in 1962, These 
investigators theorized that fluid beds could be accurately 
described by the application of the fundamental properties 
of bubbles rising in a fluidized medium. Their work and a 
recent text by Kunii and Levenspiel [15] , presents a compre- 
hensive survey of the present state of knowledge in this 
area. 

The second generation models, in general, represent the 
fluid bed as a three phase system. The dense phase, is split 
to include a cloud phase which represents solids directly 
influenced by rising bubbles. Entrapment of solids in the 
wakes of rising bubbles is also considered. Information 
on bubble size and rise velocity, cloud development, gas 
interchange rates etc. , derived from independent investiga- 
tions, are considered in model development. A representa- 
tive group of the most noteworthy second generation models, 
as compiled by Grace [16], are listed in Table I. 

The Bubbling Bed Model, developed by Kunii and Levenspiel 
[21] was shown to be semi-successful in correlating experi- 
mental reaction data and is a useful design tool. However, 
it utilizes an effective or average bubble diameter to 
describe the bubble phase. This ignores the effect of 
bubble growth on the operation of the bed and, therefore, 
seems to be an oversimplification. An improvement can be 
found in the models of Toor and Calderbank [22] , Yates et al 
[23], and Wen and co-workers [1,24], where changing bubble 
size is introduced. 
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TABLE I 



SECOND GENERATION FLUIDIZED BED MODELS 

ORCUTT et al [17] 

HOVMAND and DAVIDSON [18] 
PARTRIDGE and ROWE [19] 

CHIBA and KOBAYASHI [20] 
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Wen and Kato [24] investigated reactions using their 
Bubble Assemblege Model and found the model predicts conver- 
sions and concentration profiles reasonably well. This same 
model in modified form was used by Wen and Yoshida [1] to 
correlate results of a non-catalytic ore roasting reaction. 
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III. MODEL DEVELOPMENT 



The first step in the development of a model for a 
fluidized bed is the study of the qualitative nature of 
fluidization. In Figure 2, the three possible states of a 
bed of solids having a gas flowing through it are shown. 

At low gas flow rates, Figure 2a, the solids are in a fixed 
bed condition. The gas passes through the interstitial 
voids of the solids. As the gas flow is increased, a point 
will be reached where the pressure drop of the gas is suf- 
ficient to support the weight of the bed. At this gas flow 
rate, the bed becomes fluidized and has many of the proper- 
ties of ordinary liquids. In this state. Figure 2b, the 
bed is characterized by a minimum fluidization velocity, 
voidage, and height. 

Increasing the gas flow above the minimum fluidization 
rate causes the bed to transition to a state of vigorous and 
violent bubbling as shown in Figure 2c. This state is of 
the greatest importance for it represents the condition of 
most industrial fluid bed reactors. 

In modeling the bubbling fluidized bed, two assumptions 
are universally accepted; 1) isothermal operation and 2) the 
flow through the bubble phase represents all gas flow above 
that required for the minimum fluidization, i.e., (U-UMF) . 
These assumptions were used in the model developed in this 
investigation. This model represents a modified form of the 
Bubble Assemblege Model proposed by Wen and Yoshida [1]. 
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Figure 2. States o£ the Gas-solid System. 
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A. OVERALL BED MODEL 



The fluidized bed was divided into a number of backmixed 
reactors connected in series as shown in Figure 3. The 
height of each section corresponds to the average bubble 
diameter at that level. The model thus accounts for chang- 
ing bubble diameter within the bed. 

1 . Bubble Size and Section Height Calculation 

The bubble diameter was calculated using the rela- 
tion of Kobayashi [25] who found that experimental data on 
bubble sizes could be correlated by the equation 



DB 



0. 684ppPp 




(3) 



To account for the initial bubble size at the surface of 
the distributor Kobayashi 's relation was modified to 



DB = 0.684ppDp 




(4) 



where Do is the initial bubble size at the distributor 
plate. This quantity is difficult to estimate for beds 
having porous plate distributers and therefore must be 
estimated. For perforated plate distributors, however. Do 
can be estimated based on the work of Davidson and 
Harrison [14] . It was shown by these investigators that 
the bubbles produced in fluidized beds have essentially 
the same properties as bubbles produced in ordinary liquids 
having a small viscosity. On this basis Do can be calcu- 
lated from the equation developed for predicting the 
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Figure 3. Sectionalized Model of the Fluid Bed Reactor. 
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diameter of bubbles produced in ordinary liquids at a 
single orifice 



^ |-n5— I g 



( 5 ) 



The averaging procedure used in calculating the indi- 
vidual section heights is graphically depicted in Figure 
4. For the first section 



AHj = 



Do+XMAH^+Do 



( 6 ) 



Au _ 2 Do 
^^1 “ (2-XM) 



Repeating the procedure for the second section yields 



(7) 



AH 2 = 2Do 



(2+XM) 

(2-XM) 



In general, it can be shown the average height of the N 
section is 



( 8 ) 

th 



AHj^ = 2Do 



(2+XM) 



N-1 



(2-XM) 



N 



(9) 



The bubbles were considered to grow continuously until 
the maximum bubble diameter predicted by Harrison et al 
[26] is reached 




( 10 ) 



This maximum is attained when the upward velocity of the 
bubble equals the terminal velocity of the bed particles. 
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Figure 4. Averaging Procedure Used to Determine Section Heights. 



rH 

3 : 

<3 





24 



The calculation o£ the terminal velocity is given in the 
Appendix . 

2 . Bed Voidage and Height Calculations 

The bed conditions at minimum fluidization must be 
defined prior to the calculation of the operational bed 
characteristics. These calculations were completed using 
well defined experimental correlations that are given in 
the Appendix. 

Davidson and Harrison [14] have shown that the rise 
velocity of a crowd of bubbles in a fluid bed can be cal- 
culated from the relation 



U 



BR 



U-UMF + 0.711 




( 11 ) 



It was assumed that the expansion of the bedcover that at 
minimum fluidization is a reflection of the total volume 
of bubbles in the bed 



VB^ = (H-HMF)S^ 



( 12 ) 



and that an average bubble size can be taken as the bubble 
size at a height equal to HMF/2. The bed expansion ratio 
was calculated by combining equations (4), (11), and (12) 
to give 



H-HMF 

H 



(U-UMF) 



0.711 



gDB 

6 avg ' ^ 
12 



where 



DB = XM 
avg 



HMF 



+ Do 



(13) 



(14) 
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This equation was tested by Wen and Yoshida [1] and shown 
to correlate experimental bed heights to within 10%. 

Based on the above result, the bed voidage was calcu- 
lated as 

1 - e = (1-eMF) (15) 



B. SECTION MODEL 

The sections into which the bed was divided were treated 
as having three phases; a bubble phase, a cloud and wake 
phase and an emulsion phase. Figure 5 shows a single bubble 
and the three phases. Flow within each phase was assumed 
to be completely backmixed. This implies uniform gas con- 
centrations and uniform reactant conversion levels \>?ithin 
each phase. 

1 . Bubble Phase Description 

The equations used to describe the bubble phase 
are based on experimental results which indicate a bubbling 
fluidized bed can be treated as analogous to a bubbling 
inviscid liquid. 

The relative rise velocity of a single bubble was 
calculated from the Davies -Taylor [27] equation 



^BREL 



0.711 




(16) 



The relative rise velocity is the rise velocity of a bubble 
through a bed at minimum fluidization conditions. In a 
vigorously bubbling bed, the section of the bed ahead of 
the bubble moves upward with a velocity equal to U-UMF. 
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Figure 5. Bubble Representation. 
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The absolute rise velocity of the bed then becomes 



U 



BR 



U 



UMF + U 



BREL 



( 17 ) 



This equation was proposed by Nicklin [28] and later by 
Davidson and Harrison [14]. 

An accurate model of any section must predict the 
number of bubbles within that section. The number of 
bubbles was predicted from the bed expansion ratio according 
to the following logic 

w _ Volume of Section Comprised of Bubbles 
Volume Per Bubble 



Volume of Section X Bed Expansion Ratio 
Volume Per Bubble 



This yields 
N = 



upon 

6 



substitution 

e-cMF 

1-cMF 



The total volume of bubbles is then given by 

VB = N • DB^ 

6 



(18) 



(19) 



2 . Cloud and Wake Description 

The cloud phase was described according to the model 
of Davidson [29] which permits calculation of the cloud 
diameter according to the relation 



DC 



3 




Urrel +2UMF/CMF 

U^^^p-UMFTiMF- 



( 20 ) 



On this basis, the total volume of cloud phase in each sec- 
tion becomes 
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VC = VB 



3UMF/eMF 
^BREL - 



( 21 ) 



The cloud was assumed to have a voidage equal to the bed at 
minimum fluidization. 

The wake was described using the results of Rowe and 
Partridge [30]. These investigators have shown that the 
wake comprises approximately thirty percent of the bubble 
volume. A plot of their results is given in the Appendix 
and was used- to estimate the parameter Alpha which equals 
the ratio of the wake to bubble volume. The wake was also 
assumed to have a voidage equal to the bed at minimum 
fluidization. 

3. Emulsion Description 

The emulsion phase includes all of the remaining 
section volume 



VE = AH • S.J. - VC - VB 



( 22 ) 



having a voidage equal to the bed at minimum fluidization. 
It was assumed in this investigation that for operating 
conditions where U/UMF is large that the velocity of the 
gas through the emulsion phase is zero. This assumption 
finds substantiation in the work of Latham et al [31] and 
Kunii and Levenspiel [15] . These investigators found that 
the emulsion gas will reverse its direction (Ue=0) when 
U/UMF > 2.7 =^6.0. 

4 . Solid Movement and Gas Exchange 

The primary mechanism of solid circulation within 
the bed is by transport • in the wakes of bubbles. The solid 
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is entrained by rising bubbles and carried upward with a 
velocity equal to the rise velocity of the bubbles. During 
the bubble's movement through the bed, the solids in its 
wake are continuously exchanged with the bulk emulsion 
solids. In this investigation, the exchange coefficient 
for the solids was assumed to be infinite. 

Solids in the bulk emulsion phase move downward 
with the same volumetric flow rate as solids carried upward 
in the wakes. As a result, there is no net flow of solids 
across any horizontal plane within the bed. 

For a system in which a solid is continuously fed 
to the bottom of the bed and removed from the top, a net 
upward flow does exist and is equal to the volumetric flow 
rate of the feed. Mathematically, this model can be 
expressed as 







e ^ “^BN^BRN 
^BN 144.0 



WE 



N+1 



“^BN^BRN ^f ^ ^ . 

144.0 ■ S.J. ^^T'^BN-’ 



(23) 

(24) 



t It 

where WB(N) represents solids flowing upward to the (N+1) 

j_ -L 

section from the section and WE(N+1) represents solids 

"t It t It 

flowing downward into the section from the (N+1) 

section. 

The overall gas exchange coefficient was calculated 
using the correlation of Kobayashi et al [32]. Based on a 
unit volume of bubbles, this expression can be expressed as 

Kbe = (25) 
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C. KINETICS 



Non-catalytic gas -solid reactions in which the reacting 
particles do not change size were studied in this model. 
Examples of this type of reaction include sulfide ore roast- 
ing, reduction of iron ore, and calcination of limestone. 

The following stoichiometric equation was used to represent 
these reactions 



A(gas) + bB(solids) PRODUCTS (26) 

For this stoichiometric relation, the first order irrever- 
sible reaction rate equation can be written as 



1 dNa 

A dt 
rx 



BTT 



rx 



dNb 

~3T 



= kc Cg 



(27) 



The reaction was assumed to proceed according to the 
shrinking core reaction model [1,33]. This model assumes 
the reaction begins as the surface of the material and 
proceeds inward; the reaction zone being the surface of the 
constantly shrinking unreacted core. The particle maintains 
its structural integrity as a result of the formation of a 
completely converted inert ash which surrounds the unreacted 
core. A time history of the reacting particle is shown in 
Figure 6. It was also assumed that the overall conversion 
rate is controlled by the chemical reaction step. This 
implies that the diffusion of gaseous reactant occurs much 
faster than the reaction at the surface of the core and, 
therefore, a gas concentration gradient between the bulk gas 
phase and the core surface does not exist. 
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Figure 6. Shrinking Core Reaction Model. 



D. MATERIAL BALANCE 



1. Preliminary Assumptions 

A modified reaction rate constant based on unit 



volume of solids was defined by dividing both sides of 
Equation (27) by Vo to give 




(28) 



where 




2 



Kr 



(29) 



Kr is not a constant at a given temperature. It depends 

upon the core radius, Rc, which is a measure of the average 

level of conversion of the solids. It was assumed in this 

investigation that each of the N sections into which the 

bed was divided contains solids having a uniform conversion 

level or stated in terms of the shrinking core model having 

an equal unreacted core size. Furthermore, because of 

backmixing within each section, the solids swept upward 

from the section to the (N+1)^^ section in the wake 

of bubbles were assumed to have a conversion level charac- 

t h 

teristic of the section. This is also true for the 

downward flowing emulsion solids. These assumptions are 

shown in Figure 7. It should be noted that no real physical 

significance can be attached to the concept of having an 

average core radius within each section. It is obvious 

that the radius of the solids entering the section from 

t h 

the (N+1)^ section does not increase to that of the material 
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Figure 7. Core Size Distribution. 
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section . 



in the N 



th 



This assumption was made merely to 



simplify the mathematical representation in terms of the 

X 

section material balances. 

To further simplify the physical representation, the 
shrinking core model was used to define an average solid 
reactant concentration for any particular solid conversion 
level. This approach can be explained with the aide of 
Figure 8. It was assumed that within each particle, reac- 
tant is uniformly dispersed through the entire volume of 
the particle. Since the shrinking core model states that 
reaction takes place only at the surface of the unreacted 
core, using the above assumption it can be stated that the 
unreacted core of a partially converted particle has a 
reactant concentration equal to the initial reactant con- 
centration. The average concentration for the entire 
particle was then estimated as 



2 . Section Material Balances 

Steady state mass balances on the gaseous reactant 
were used to characterize the flow streams of the system. 
Figure 9 shows an individual se^ction and the associated 
gas flow streams. The material balances for the bubble 
phase and emulsion phase can be written as 



CS V 
o core 



(30) 



V 

o 





( 31 ) 
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Figure 8. Radial Solid Concentration Profile. 
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Figure 9. Section Gas Flow Representation. 
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CGB(N-l) 



\ 



(Kbe)^[CGB^-CGEj^] - YeNKr^CGEj^ = 0 



(32) 



Included in the reaction term are y and y which define 

c e 

the distribution of solids between the two phases of the 
section. These terms are defined as 



The material balances were derived on the basis of zero gas 
flow in the emulsion and complete backmixing within each 
phase . 

In the top reactor section, the cloud phase volume must 
be modified to account for extension of the bubbles above 
the bed surface. The bubble was described as shown in 
Figure 10. In this modification, only the fractional volume 
of the cloud phase associated with that part of the bubble 
below the bed surface is considered. The equations for 
cloud phase then become 



= Volume of Solids Dispersed in Clouds and Wakes 



Y^ = tVCj^+ aVBj^] • [l-cMF]/1728.0 



(33) 



Yg = Volume of Solids Dispersed in The Emulsion 



Ye = VEj^ • [l-cMF]/1728.0 



(34) 




where 



DIFH < 0.5 AH. 



N 



Vpj^ = 1.0472 [DIFH] ^ tS.OR^jj^- DIFH] 



(35) 
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Figure 10. Top Section Modification. 
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VBpj^ = 1.0472 [DIFH-R^^+Rg]^[2.0Rg-DIFH+R^^] (36) 

DIFH >0.5 

VBpR = V.pp- 1.0472 [AH-DIFH]^ [3. OR^^-AH+DIFH] (37) 

VBpj^=VBp- 1.0472 [AH-DIFH-R^^ + Rp] ^ [ 2 . ORg- AH+DIFH + R^^ ] 

(38) 

An overall material balance on both the solid and gas 
streams was derived from Figure 11 which shows the bed 
representation from the (N-1)^^ section down to the feed 
section. An overall balance yields 

‘>«vol lCG„-CGB^.p = W^CS^ * WE^CS^ 

^^N-l^^N-1 (39) 

3 . Solution of the Material Balance Equations 

A trial and error procedure was used fo" obtain solu- 
tions to the material balances of the system. The procedure 
can be characterized by the three major steps involved 

Step 1: Assume the exit gas stream concentration CGB(N) 

Step 2: Solve each of the section material balances to 

obtain the required initial gas concentration 
for the assumption of Step 1. 

Step 3: Compare the calculated initial gas concentra- 

tion with the known value (CGq) and continue 
the calculation until they are equal by 
adjusting the value of CGB(N). 

A detailed description of the completion of Step 2 is given 

in the next three sections. 

a. Calculation for the Top Section 

The bed with the section isolated is shown 
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Figure 11. Flow Streams for Overall Material Balance. 
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in Figure 12. Four material balances can be derived from 



this figure: 

- V<.f,Kr^CGB^ - 0 (40) 

0'be)N''®N'C«»N-CGEN' ' ■^eN*'''NCGEn = 0 (41) 

' *'f“o " ■ ^N-l^N-1 f''®) 

bUvol'‘^®®N-r‘^®®N’ ' ™N-i“n- 1 ■ 'V®N ' f'*®’ 



Equations C40) and (41) are the gas reactant material 

balances for the bubble and emulsion phase, Equation (42) 

th 

is the overall material balance over the (N-1) section 

to the feed section, and Equation (43) is the balance for 

t h * 

both streams in the section. 

The solution of these equations requires knowledge of 

i- 

the average core size of material in the section for 

evaluation of Kr, the reaction rate constant. This is not 
known and therefore a second trial and error solution pro- 
ceeds as follows . 

Step 1: Guess Rcore and calculate Kr from Equation (29) 

and Vcore. 

Step 2: Calculate CGE(N) and CGB(N-l) from Equations 

(41) and (40) . 

Step 3: Calculate CS(N) from Equation (30). 

Step 4: Calculate CS(N-l) from Equation (42). 

Step 5: Repeat the calculations until the section 

material balance Equation (43) is satisfied. 
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Figure 12. Top Section Material Balances. 
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When Step 5 has been completed, the solution proceeds to 
the (N-1)^^ section. 

b. Calculations for the (N-1)^^ Section 

The middle bed sections are characterized by 
the three standard material balance equations which can 
be written as 

(Kbe)^.^VB^.^[CGB^.j-CGE^.jl - YeN-lKl-N-iCGE^-rO (45) 



bU , [CG -CGB„ ,1 = W.CS + WE., , CS,, , 
vol^ o N-1 f o N-1 N-1 



+ WBj^_2CSj^_2 0 



(46) 



Solution of these equations does not require a trial and 
error procedure sin.ce the average core size of material in 
the (N-1)^^ section is fixed by the results of the calcu- 
lations on the (N)^^ section through 

CS Vcore, 



CS 



'N-1 



N-1 



V, 



(47) 



This equation defines Rcore and Kr for the (N-1)^^ section. 
Solution then proceeds as follows 



Step 5: CGE(N-l), CGB(N-2) are calculated from 

Equations (45) and (44). 

Step 6: CS(N-2) is calculated from Equation (46). 

This procedure is repeated until the feed section is reached. 
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c. Calculation for the Feed Section 

The material balances around the feed section 
can be written as 

^vol^^^^o'^^^l^ CKbe)iVBi[CGE^-CGB^] 

- Y^^Kr^CGB^ = 0 (48) 

(Kbe)iVBi[CGBi-CGEi] - Ye^Kr^CGE^ - 0 (49) 

As is the case for the middle sections, a trial and error 
procedure is not required for their solution on this basis, 
solution was completed as follows 

Step 7: Rcore is calculated from CS(1) obtained from 

the calculations on section two. 

Step 8: Kr and CGE(l) are calculated from Equations 

(29) and (49). 

Step 9: CG^ is calculated from Equation (48). 

The initial gas concentration calculated in Step 9 is then 
compared with the known value of CG^. This provides a 
test of the validity of the assumed value of the exit gas 
concentration CGB(N). If the calculated CG^ does not equal 
the experimental value, CGB(N) is adjusted and the entire 
calculation repeated. When the equality exists, the material 
balances for each section and the overall material balance 
for the bed have been satisfied and the solution completed. 
The computer flow chart for the solution of the material 
balances is given in Figure 13. 
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Figure 13. Flow Chart for the Solution 
of the Material Balances. 
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E. THE COMPUTER MODEL 



The program developed to model the fluid bed reactor is 
given in the COMPUTER PROGRAM SECTION. The program was 
designed to be completely self-contained, that is, no 
external calculations are required. The input variables 
are listed in Table II. They represent common design 
variables and therefore would normally be well defined. 

In Table III the program is divided into its three major 
sections and the calculations completed in each are listed. 



\ 
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TABLE II 



COMPUTER INPUT VARIABLES 


1. 


Reactor Diameter 




2. 


2 

Number of Holes per in. in 
Distributor (0 for Porous 


the 

Plate) 


3. 


Weight of Solids in the Bed 




4. 


Diameter of the Bed Solids 




5. 


Radius of Reactant 




6. 


Density of Bed Solids 




7. 


Density of Reactant 




8. 


Bed Temperature 




9. 


Stoichiometric Constant 




10. 


Kinetic Rate Constant 




12. 


Molecular Weight of Solid Reactant 


13. 


Solid Feed Rate 




14. 


Initial Solid Concentration 


' 


15. 


Initial Gas Concentration 




16. 


Superficial Gas Velocity at 
Temperature 


Room 


17. 


Superficial Gas Velocity at 
Temperature 


the Bed 


18. 


Alpha (Volume of Wake/Volume 


Bubble) 


19. 


Initial Section Height 




20. 


Initial Guess of CGB(N) 




21. 


Increment of Gas Search 
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MAJOR COMPUTER PROGRAM SECTIONS 
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Solid Flow Streams 



IV. PRESENTATION AND DISCUSSION OF RESULTS 



The model developed in the previous section was used to 
study a sulfide ore roasting reaction following the stoichio- 
metric equation 

2ZnS + 3 O 2 2AnO + 2 SO 2 (50) 

The experimental results obtained by Yagi et al [34] were 
used for evaluating the model. The performance of the model 
was judged by its ability to predict experimental conver- 
sions. This presentation is divided into two sections. In 
the first, the overall results are presented predicted by 
the model of Wen and Yoshida [1] . In the second section, 
the effects of the major assumptions used in the development 
of the model are investigated. This section represents a 
detailed study of the model performance. 

A. OVERALL RESULTS: GENERAL MODEL PERFORMANCE 

In Table IV, the experimental data of Yagi et al [34], 
the calculated conversions obtained by Wen and Yoshida [1] , 
and those from the model developed in this investigation 
are tabulated. In these runs, the height of the first 
section AH^ was assumed to be 1.0 cm. This assumption was 
proposed by Wen and Yoshida for reactors having porous plate 
distributors. The validity of this assumption is discussed 
in another part of this section. The input data for these 
runs are given in Table VIII of the Appendix. 
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Run 

No. 

1 

2 

3 

4 

6 

7 

10 

11 

12 

13 



TABLE IV 



COMPARISON OF OVERALL RESULTS WITH LITERATURE DATA 



Conv. 
Exp . 


Conv. 

Model 


% 

Error 


Conv. 

WEN 


% 

Error 


99.4 


99.9 


0.59 


99.8 


0.43 


97.2 . 


99.9 


2.77 


98.5 


1.34 


88.7 


99.4 


10.75 


90.0 


1.47 


91.0 


99.2 


8.31 


93.0 


2.18 


86.5 


97.2 


11.05 


86.9 


0.46 


93.5 


99.7 


6.22 


95.0 


1.61 


85.4 


97.1 


12.01 


- 


- 


80.6 


91.2 


11.62 


- 


- 


72.4 


81.8 


11.49 


- 


- 


91.7 


99.9 


8.17 


- 


- 
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The data shows that the model predicts conversions which 
are significantly higher than the experimental values. For 
experimental runs having relatively high conversions (>90%), 
the material is essentially completely reacted within the 
bed (>99% conversion) . The average percent error of the 
calculated conversions for all runs was 8.3%. The maximum 
error was 12.01% for Run 10 and the minimum was 0.59% in 
Run 1. The model did not predict for any of the runs, a 
conversion less than the experimental value. Low percentage 
error for runs having higher conversion can therefore be 
expected. 

The calculated conversions reported by Wen and Yoshida 
show much better agreement with the experimental data. The 
average percentage error for their runs was 1.25%. These 
investigators did not report results for experimental runs 
in which low conversions were found and for which the model 
proposed in this investigation had the highest error. On 
the basis of the six runs reported by Wen and Yoshida, the 
error for the proposed model becomes 6.60% . 

B. EFFECTS OF THE MAJOR ASSUMPTIONS 

1 . Modified Section Height Calculation 

The height of each section was originally calculated 
by taking the average bubble diameter at a particular height 
as described in the MODEL DEVELOPMENT SECTION. When this 
is done, a portion of the cloud surrounding the bubble 
actually extends beyond the section boundaries as shown in 
Figure 14. The volume of the cloud phase associated with 
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Figure 14. Section Height Based on Bubble Size. 



Cloud 




Hatch areas show volume of cloud 
phase which extends above and 
below section boundaries 
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each section, therefore, is slightly larger than the correct 
value. Furthermore, this portion of cloud volume is counted 
twice when succeeding sections are considered. 

Figure 15 shows a plot of the gas concentrations for 
the bubble and emulsion phase as a function of height above 
the distributor. It can be seen that the bubble gas concen- 
tration is considerably greater than the concentration of 
gas in the emulsion phase. As a result, solids in the 
bubble wake and cloud react with gas of higher concentration 
than emulsion solids. It was theorized that the inaccura- 
cies in the prediction of the volume of the cloud phase 
resulting from the section height calculation, might con- 
tribute to the prediction of a larger overall conversion. 

To test this theory, the section height calculation was 
modified to include the entire cloud volume. The equations 
necessary to make this modification and their derivation 
are given in the Appendix. 

The conversion results obtained from the modified pro- 
gram are given in Table V. No significant change in the 
predicted conversions was found. The overall percent error 
was 8.30% which equals that found for the original program. 
The average error for the six runs reported by Wen and 
Yoshida was found to be only slightly higher at 6.64%. The 
maximum change in predicted conversion occurred for Run 12 
where the conversion calculated from the modified program 
was found to be 1.6% higher than the original calculation. 
These results show that the error in cloud phase volume 
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Figure 15. Gas Concentration Profiles for Run 12; 
AH^ Based on Cloud Size. 
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TABLE V 



PREDICTED CONVERSIONS FOR MODEL WITH SECTION HEIGHT 
BASED ON THE CLOUD DATA 



Run 

No. 




Conv. 

Exp. 




Conversion 


% 

Error 

B 


A 


B 


1 




99.4 




99.9 


99.9 


0.59 


2 




97.2 




99.9 


99.9 


2.77 


3 




88.7 




99.4 


99.4 


10.73 


4 




91.0 




99.2 


99.3 


8.36 


6 




86.5 




97.2 


97.6 


11.10 


7 




93.5 




99.7 


99.7 


6.24 


10 




85.4 




97.1 


97.2 


12.15 


11 




80.6 




91.2 


91.9 


12.29 


12 




72.4 




81.8 


83.1 


10.67 


13 




91.7 




99.9 


99.8 


8.17 


A: 


Section 


Height 


Based 


on Bubble 


Diameter 




B: 


Section 


Height 


Based 


on Cloud 


Diameter 
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mmmm 












resulting from calculating the section heights on the basis 
of average bubble size is not a significant factor. In the 
rest of this study, however, the modified section height 
calculation was retained. 

2 . Backmixing Considerations 

In the development of this model, it was assumed 
that the solids within a particular section have an equal 
average concentration and conversion level. This concept 
implicitly assumes that every particle fed to the reactor 
has a residence time which is equal to that of any other 
particle. The mean residence time being defined as 



t = 




(51) 



Investigators [15,35] studying the flow patterns and 
residence times in fluidi zed ' beds have shown that a residence 
time distribution (RTD) actually exists for the feed 
particles. This concept implies that some feed solids have 
a residence time which is far less than the average, while 
others have residence times greater than the average. As 
a result, the product of a reaction occurring in a fluidized 
bed does not have a uniform conversion level. Kunii and 
Levenspiel [15] proposed that the RTD can be expressed by 
defining an exit age distribution function E(t) as 

E(t) = (52) 

t 

where E(t)dt is the fraction of solids staying in the bed 
for a time between t and t+dt. This expression has been 
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shown to correlate experimental residence times very 
accurately . 

The normalized exit age distribution function is plotted 
in Figure 16. The shaded area represents the fraction of 
material having a residence time equal to t+10|. This 
fraction was found to be only 1 . 1 % of the total bed 
material. This observation indicates that the assumption 
that all the solids in the bed have a residence time equal 
to t” might be a gross oversimplification leading to poor 
calculated results. 

To test this theory, overall bed conversions based on 
a non-uniform bed residence time were calculated using 
average gas concentrations predicted by the model. The 
equations used were those given by Kunii and Levenspiel 
jl^ . According to these authors, the mean conversion of 
the product stream can be defined as 



XB = 3 


I M - 6 


t 


2 

+ 6 


t 


3 


1.0 - 




^ ) 


T 




T 







(53) 



where the time for complete conversion is equal to 



T- = P'PDP (54) 

2bkcCg 

The overall average gas concentration ^ was calculated as 
a weighted average of the section gas concentrations. 

N 

Ug = .L__ (55) 

t 

The results of the calculation are presented in Table VI. 
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Figure 16. Exit Age Distribution. 
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TABLE VI 



RESULTS OF CONVERSION CALCULATIONS 
FROM BACKMIXED EQUATIONS 



Run 

No. 


Conv. 
Exp . 


Conv. 

WEN 


Conv. 

Model 

B.M. 


% 

Error 
for B.M. 
Convers . 


1 


99.4 


99.8 


99.5 


0.1 


2 


97.2 


98.5 


98.9 


1.75 


3 


88.7 


90.0 


95.6 


7.78 


4 


91.0 


93.0 


94.9 


1 

4.29 


6 


86.5 


86.9 


87.8 


1.60 


7 


93.5 


95.0 


97.1 


3.85 


10 


85.4 


- 


1 

8 7.9 


2.93 


11 


80.6 


- 


75.7 


6.08 


12 


72.4 


- 


61.2 


15.47 


13 


91.7 


- 


98.2 


7.09 
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The conversions calculated from the backmixed equations 
show considerably better agreement with experimental results. 
The overall average error was reduced to 5.081 and the error 
for the first six runs to 3.211. These compare with errors 
of 8.301 and 6.64% obtained from the model. These results 
indicate that neglecting the existence of a residence time 
distribution for the bed particles causes the prediction of 
conversions which are high when compared with experimental 
results . 

The Bubble Assemblage Model proposed by Wen and Yoshida 
does account for the non-uniform residence times of the bed 
particles. The close agreement of their results with experi- 
mental data is a further indication that the assumption of 
uniform solid conversion levels in each section incorporated 
in the model proposed in this investigation is an over- 
simplification. 

3 . Analysis of the Initial Reactor Section 

One of the major assumptions required in the analysis 
of the model performance was the assumption that the height 
of the first section is equal to 1.0 cm. This assumption is 
necessary because the reactor used by Yagi had a porous plate 
distributor for which an accurate calculation of initial 
bubble size is not possible. The effects of this assumption 
on their model's performance were tested by Wen and Yoshida. 
They concluded that the calculated conversions were not sen- 
sitive to the height of the first section. Gas and solid 
reactant concentration profiles given by these investigators, 
however, show very large changes occurring in the area 
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directly above the distributor. This observation indicates 
that the conclusion that the initial section height does 
not effect the overall results might not be adequate for 
all reaction conditions. In this section, the results of 
a detailed study of the first reactor section are given. 

It was found that essentially all the reaction occurs 
in the first section of the bed. Data supporting this 

/ 

observation are given in Table VII. These data show that 
the concentration of solids leaving the first section is 
reduced to very nearly the final product concentration 
level. This observation is not surprising. Many investi- 
gators [36,37] have reported that the primary reaction zone 
is the area directly above the distributor. This fact 
emphasizes the necessity to accurately model the initial 
section . 

Wen and Yoshida have suggested that the only criterion 
required in determining the effects of the first section 
height on the overall model performance is the kinetic 
speed of the reaction, i.e., the value of the rate constant 
kc. They conclude from their results that the overall con- 
version for a fast reaction is not effected by the value of 
and, therefore, no effect will be seen for reactions 
having low rate constants. 

In Figures 17 and 18 gas concentration profiles for 
Runs 1 and 12 are plotted. These two runs represent the 
extremes in experimental conversions which were investigated. 
Run 1 is kinetically the faster of the two; having a kc = 

6.56 X 10 ^ ft/sec compared v\rith kc = 0.328 x 10 ^ ft/sec 
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TABLE VII 



SOLID CONCENTRATION DATA 



Run 

No. 


Initial 
Cone. , 
mole/£t . 


Avg. Cone. 

Leaving 
l^t Seetion 


Avg. Cone. 

Leaving 
Top Seetion 


1 


'2.22 


4.3 X 10‘^ 


2.5 X lO""^ 


2 


2.22 


1.0 X 10"^ 


7.6 X lO""^ 


6 


3.58 


8.7 X 10"^ 


8.6 X 10'^ 


10 


3.58 


1.0 X 10"^ 


9.9 X 10"^ 
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Figure 17. Gas Concentration Profiles for Run 1; 
AH^ Based Cloud Size. 
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Figure 18. Gas Concentration Profiles for Run 12; 
AH^ Based on Cloud Size. 
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for Run 12. The figures show that gas concentrations for 
the kinetically slower reaction change more rapidly in the 
initial height of the bed than do the concentrations for 
the more rapid reaction. The emulsion gas concentration 
showed a very rapid drop for Run 12. At a height of 0.5 in. 
from the distributor, a 711 decrease in emulsion gas^ concen- 
tration was found. This compares to a 27% change in the 
emulsion gas concentration for Run 1. 

The effect on the overall conversion resulting from a 
change in the initial section height for Runs 1 and 12 is 
shown in Figure 19. As suggested by the gas concentration 
profiles, Run 12 shows a sensitivity to AH^ , while Run 1 
does not. This trend is not expected if Wen's and Yoshida's 
conclusion that slow reactions are not sensitive to the value 
of AH^ is accepted. 

In Run 12, the reactant feed rate and the initial reac- 
tant concentration of the feed were substantially larger 
than Run 1. The fact that the conversions in Run 12 were 
sensitive to the value of AH^^ suggests that the feed rate 
and initial concentration as well as the value of the rate 
constant are important parameters in determining the effect 
of AHj^ on the overall bed conversion. 

Runs 10, 11, and 12 were used to study this effect. 

These runs represent a series in which all parameters 
except the feed rate are constant. The feed rate increases 
for this series of runs. 

In Figure 20, normalized average gas concentrations are 
plotted for the three runs. As the feed rate of the solid 



66 



95 

90 

85 

80 

75 

70 

65 



0 — 0 



0 



0 



Run 1 




:ure 19. 



Conversion Versus Runs, 1, 12. 



67 



Normalized Average Gas Concentration 
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Figure 20. Normalized Gas Concentration 
Profiles Runs 10, 11, 12. 
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to the bed increases, conversion and average gas concentra- 
tions decrease as would be expected. All the runs show 
rapid gas concentration changes near the distributor plate 
which suggests that the value of AH^ will effect the overall 
calculated conversion. This was found to be true as shown 
in Figure 21. Run 12 showed the greatest sensitivity to 
a change in the initial section size. The conversion for 
Run 12 decreased by 21% over the range of . For Run 11, 
the change was only 1S% and for Run 10, it was only 5%. 

These results show that the feed rate does effect the sensi- 
tivity of the model to the initial value of AHj^ . 

The complete meaning of this conclusion is not clear, 
but some speculation concerning the reasons for this feed 
rate dependency can be made. The fact that conversion was 
found to be dependent on the. value, of AH^ can be interpreted 
to mean that the model will predict experimental conversions 
accurately if correct values of AH^ are known. This logic 
ignores the fact that the model might still predict incorrect 
conversions due to its not accounting for a residence time 
distribution of solids within the bed, an effect which was 
previously shown to be significant. 

From Figure 21, the values of AH^ , for which the model 
will predict the correct experimental conversions can be 
read. The values are: Run 10, > 1.0; Run 11, 0.88; and 

Run 12, 0.83 in. The interesting fact is that AHj^ is 
different for the three runs and decreases as the feed flow 
rate is increased. One would expect that for the same gas 
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Figure 21. Conversion Versus AH^ Runs 10, 11, 12. 
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flow conditions, bubbles of equal size would form at the 
distributor. The trend, however, suggests that different 
size bubbles are being formed at the distributor and that 
their size is dependent upon the reaction conditions in 
the bed. 

Recent investigation [37] into the mechanics of the 
bubble formation at distributors in fluidized beds have 
shown that very small bubbles are initially formed, which 
because of their large inertia swirl in turbulent jets just 
above the distributor surface. These small bubbles rapidly 
coalesce and form large, stable bubbles at the tips of the 
gas jets. If this model of bubble formation is correct, a 
dependence between the size of the bubbles formed at the 
jets and the reaction conditions is not surprising. 

Consider the effects of increasing the feed flow rate 
to the bed. This action increases the tendency for reaction 
and thus the disappearance of gas from the bubble phase. 

As a result, two effects will occur: 1) smaller, stable 

t 

bubbles will be formed by coalesence and, 2) sharper gas 
concentration drops will occur. These trends are both 
present in the data plotted in Figures 20 and 21; a feet 
which supports the speculative logic used in the analysis. 

The conclusion that can be drawn from these results is that 
the uncertainty in the existing mechanism of bubble forma- 
tion could possibly be the controlling factor in the develop- 
ment of a successful mathematical model of fluidized bed 
operations . 
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V. CONCLUSIONS AND RECOMMENDATIONS 



The major conclusions of this investigation are: 

1. The conversions predicted by the proposed model are 
significantly higher than the experimental values of Yagi 
et al [34] . The average percent error of the predicted 
conversions was 8.31. 

2. The results are not dependent on the method used for 
calculating the section heights. Conversions obtained from 
the model when the section heights were based on the cloud 
diameter were essentially identical to the conversions 
predicted by the model when the section heights were based 
on the bubble diameter alone. 

3. The assumption that the bed solids have a uniform resi- 
dence time equal to t", which ignores a solid residence time 
distribution, contributes to the calculation of conversions 
which are higher than experimental values. Calculations 
using the residence time distribution concept and the 
average gas concentration predicted by the model show a 
significant improvement in the overall results. For these 
calculations, the average percentage error was reduced to 
5.081. Modification of the model material balances to 

include the residence time distribution concept is recommended. 

■ 0 - 

4. Investigation into the effects of having to assume the 
height of the initial reactor section indicate that the 



72 



I 



I 



i 



i 




overall conversion results are dependent upon the assumed 
value (AH^) . Large changes in the solid and gas concen- 
trations in the initial section support this conclusion. 

The sensitivity of the predicted conversions to the value 
of AH^ was found to be a function of the kinetic rate con- 
stant and the feed rate of solid reactant. The trends 
indicate that the sensitivity increases for decreasing rate 
constants and for increasing feed flow rates. 

The model proposed in this investigation represents 
the first step in the development of a mathematical model 
for use in the design and study of fluidized bed reactors 
for the combustion of solid wastes. Modification to 
include the concept of a solid residence time distribution 
is required. Furthermore, complete understanding of the 
mechanism of bubble formation at the surface of the dis- 
tributor in a fluidized bed is essential in light of the 
observed dependence between the overall conversion and the 
initial section height. For use in the study of combustion 
type reactions, modification for handling the kinetics of 
shrinking particles must also be made. In this regard, the 
concepts proposed by Kunii and coworkers [15,38] are 
recommended. 
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APPENDIX A 



TERMINAL VELOCITY CALCULATION 



The terminal velocity of the bed particles was calcu- 
lated using the equations derived by Leva [39] . 

2 

(pp-pp)gDp 

■■ 



18 P 



Re < 2.0 



(56) 



. 1 ron„1.14„.714. „ ^ ..714 
0.152Dp g (PP-P„) 

n = 2 

^T ,.428 .285 

P P„ 



20<Re<500 (57) 



= 



’3gDp(pp-Pg) 



n. 5 



Re > 500 



(58) 
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APPENDIX B 



CHARACTERISTICS AT MINIMUM FLUIDIZATION 



The bed voidage at minimum fluidization was calculated 
from the data of Agarwal and Storrow [41] . The cMF was 
found by these investigators to be a function of the 
particle diameter and the bed material. Their data for soft 
brick particles is shown in Figure 22. This data was fitted 
to three straight lines. The correlation coefficients of 



The equations are 




Dp £ 


0.003 




cMF = 


-38.0Dp + 0.613 


(59) 


0.003 


< Dp < 0.006 




cMF = 


015.8 Dp + 0.542 


(60) 


Dp £ 


0.006 




= -4. 


18 Dp + 0.474 


(61) 



The minimum fluidization velocity was calculated from 
an equation given by Kunii and Levenspiel [15] . For small 
particles this equation is 



UMF 



2 

Dp (pp-pg)g 

1650y* 



(62) 



The gas viscosity and density required for this calculation 
were calculated at the bed temperature and atmospheric 
pressure from the equations 

PM 

pg = ^ IDEAL GAS LAW (63) 
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Figure 22. EMF Data. 
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Ti + 1.47 T 

yii = yi2|TT 



B 



+ 1.47 Tg 



(64) 



mix 



z y.(M.)- 



(65) 



Equations (64) and (65) were obtained from Perry’s Handbook 
for Chemical Engineers [41] . 

The height of the bed at minimum fluidization was cal- 
culated from the EMF by the equation 

12 H_ 



HMF = 



1.0-eMF 



( 66 ) 



where 



144 W, 



= c 

o p S 



(67) 



P t 
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Vw/(Vw+V ) = a/(l+a) 



APPENDIX C 




dp (mm . ) 



Figure 23. Plot of Data for Determination 
of the Wake Volume. 
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APPENDIX D 



TABLE VIII 
INPUT DATA 



Run 

No. 


‘’P -3 
X 10 
in. 


Temp . 
op 


X 10 ^ 
Ib/sec 


U 

ft/sec 


kc 

X 10’^ 

ft/sec 


CGo 

X 10 ^ - 
mole/f 


CSo 

mole/£t^ 


1 


4.33 


1652 . 


6.65 


.6562 


6.56 


1.361 


2.22 


2 


4.33 


1652. 


9.941 


.6562 


6.56 


1.361 


2.22 


3 


3.50 


1292. 


3.674 


.6365 


. 328 


1.648 


3.577 


4 


3.50 


1382. 


7.458 


.6365 


.984 


1.561 


3.577 


6 


3.50 


1382. 


11.023 


.6365 


.984 


1.561 


3.577 


7 


3.50 


1472. 


7.349 


.6365 


1.97 


1.486 


3.577 


10 


3.50 


1292. 


7.349 


.6365 


.3281 


1.648 


3.577 


11 


3.50 


1292. 


11.023 


.6365 


.3281 


1.648 


3.577 


12 


3.50 


1292 . 


14.680 


.6365 


.3281 


1.648 


3.577 


13 


3.50 


1382 . 


3.895 


.6693 


.984 


1.561 


3.577 
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APPENDIX E 



MODIFIED SECTION HEIGHT CALCULATION 

The cloud diameter surrounding a bubble of size DB can 
be calculated from the relation proposed by Davidson [29] 



DC = DB 



+ 2 UMF/eMF 
- UMF/eMF 



nl/3 



( 68 ) 



In order to define a section height (AH) on the basis of 
the cloud diameter, an averaging procedure similar to that 
employed previously was used. Using this procedure, the 
average section bubble diameter was defined as 



N 

XM E H + Do 

DB^ = (69) 

^ 1 - XM/2 



On this basis, the section height based on the cloud 
diameter can be written as 



When 

each 



"»N = 



+ 2 UMF/eMF 
Ugj^ - UMF/eMF 



nl/3 



(70) 



this formulation is used, the number of bubbles in 
section becomes 



6 (H-HMF) AH., 

N = i ^ 

ttH DBj^ 



(71) 
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